
library(haven)
library(stargazer)
library(foreign)
library(stringr)
library(data.table)
library(gdata)
library(ggplot2)
library(lfe)
library(reshape2)
library(Rmisc)
library(xtable)
library(multcomp)
library(dplyr)
library(cowplot)
library(rdd)
library(bit64)
library(DataCombine)
library(rdrobust)
library(zoo)
library(RColorBrewer)
rm(list=ls())



specify_decimal <- function(x, k) format(as.numeric(round(x, k), nsmall=k))
ihs <- function(x) log(x + sqrt(x^2+1))

mod_stargazer <- function(est) {
  capture.output(est)
}


### Load Municipal Level Data
### Linked consituencies
load(file="vrn_b.Rda")
load(file="vrn_b_year0.Rda")
### County data
load(file="exp_full.Rda")
### Candidate data
load(file="candall.Rda")



#########################################################
#############     MAIN TEXT           ###################
#########################################################

###########################################
########     TABLE  1      ################
###########################################
exp_full_2015<-subset(exp_full, year==2015)

exp_full_2015<-subset(exp_full_2015, unit_type!="внутригородская территория города федерального значения" &  unit_type!="")

cand_unit<-candall[,list(num_candidates=uniqueN(vibid),
                         percentur=max(ur[wonofficial==1],na.rm=TRUE)
                         
),by=c("vrn","unit_type")]
cand_unit<-cand_unit[,list(num_candidates=mean(num_candidates,na.rm=TRUE),
                           percentur=mean(percentur)*100
                           
),by=c("unit_type")]


cand_unit<-subset(cand_unit, unit_type!="" & is.na(unit_type)==FALSE)
exp_full_2015_col<-exp_full_2015[,list(N=.N,
                                       totalexpend=mean(totalexpend,na.rm=TRUE)/1000,
                                       totalincome=mean(totalincome,na.rm=TRUE)/1000,
                                       transfers=mean(transfers,na.rm=TRUE)/mean(totalincome,na.rm=TRUE)*100,
                                       allpopulation=mean(allpopulation,na.rm=TRUE)/1000,
                                       totalterritory=mean(totalterritory,na.rm=TRUE)/1000,
                                       total_contracts=mean(total_contracts,na.rm=TRUE)/1000000,
                                       total_construction_contracts=mean(total_construction_contracts,na.rm=TRUE)/1000000                     
),by=c("unit_type")]
exp_full_2015_col<-merge(exp_full_2015_col,cand_unit,by=c("unit_type"),all.x=TRUE,all.y=FALSE)

exp_full_2015_col$unit_type=NULL
exp_full_2015_col<-t(exp_full_2015_col)

exp_full_2015_col<-cbind(exp_full_2015_col[,2],exp_full_2015_col[,3],exp_full_2015_col[,1],exp_full_2015_col[,4])

row.names(exp_full_2015_col)<-c("N","Total Expenditures (mil. rub)","Total Revenue (mil. rub)","Dependence on Subsidies (%)","Population (ths.)","Territory (ths. hectars)","Mayoral Procurement (mil. rub)","Mayoral Construction Procurement (mil. rub)","No. Candidates per Election","Won by UR Candidate (%)")


exp_full_2015_col<-as.data.frame(exp_full_2015_col)

###########################################
########     TABLE  2      ################
###########################################

SummaryTables <- data.frame(title=numeric(0),allcands= numeric(0),businesscands= numeric(0),nonbusinesscands=numeric(0))

SummaryTables[1 ,] <- c("(1) Total No.", 
                        uniqueN(candall$vibid),
                        uniqueN(candall$vibid[candall$businessperson==1]),
                        uniqueN(candall$vibid[candall$businessperson==0]))

SummaryTables[2 ,] <- c("(2) Mean Age", 
                        mean(candall$age,na.rm=TRUE),
                        mean(candall$age[candall$businessperson==1],na.rm=TRUE),
                        mean(candall$age[candall$businessperson==0],na.rm=TRUE))                  
SummaryTables[3 ,] <- c("(3) Female (%)", 
                        mean(candall$female,na.rm=TRUE)*100,
                        mean(candall$female[candall$businessperson==1],na.rm=TRUE)*100,
                        mean(candall$female[candall$businessperson==0],na.rm=TRUE)*100)        

SummaryTables[4 ,] <- c("(4) Education Level", 
                        mean(candall$edu_level,na.rm=TRUE),
                        mean(candall$edu_level[candall$businessperson==1],na.rm=TRUE),
                        mean(candall$edu_level[candall$businessperson==0],na.rm=TRUE)) 

SummaryTables[5 ,] <- c("(5) Political Independents (%)", 
                        mean(candall$ind,na.rm=TRUE)*100,
                        mean(candall$ind[candall$businessperson==1],na.rm=TRUE)*100,
                        mean(candall$ind[candall$businessperson==0],na.rm=TRUE)*100)               

SummaryTables[6 ,] <- c("(6) Ruling Party Member (%)", 
                        mean(candall$ur,na.rm=TRUE)*100,
                        mean(candall$ur[candall$businessperson==1],na.rm=TRUE)*100,
                        mean(candall$ur[candall$businessperson==0],na.rm=TRUE)*100) 

SummaryTables[7 ,] <- c("(7) Elections Won (%)", 
                        mean(candall$wonofficial,na.rm=TRUE)*100,
                        mean(candall$wonofficial[candall$businessperson==1],na.rm=TRUE)*100,
                        mean(candall$wonofficial[candall$businessperson==0],na.rm=TRUE)*100) 

SummaryTables$allcands<-prettyNum(specify_decimal(as.numeric(SummaryTables$allcands),1),big.mark=",")
SummaryTables$businesscands<-prettyNum(specify_decimal(as.numeric(SummaryTables$businesscands),1),big.mark=",")
SummaryTables$nonbusinesscands<-prettyNum(specify_decimal(as.numeric(SummaryTables$nonbusinesscands),1),big.mark=",")
SummaryTables[1,]<-gsub("\\.0","",SummaryTables[1,])

colnames(SummaryTables) <- c(" ","All Cands.","Businesspeople","Non-Businesspeople")

###########################################
########     TABLE  3      ################
###########################################

####### Panel A

est1<-felm(natecon_avg ~ business_win+ natecon_year0 + totalexpend_log_year0 |  factor(unit_type) | 0 | region, data=vrn_b, exactDOF=TRUE)

est2<-felm(natecon_avg ~ business_win + natecon_year0 + totalexpend_log_year0  + transfers_level_year0 + turnout + voterlist_log + numcandsperelection_log + incumbent |  factor(unit_type) + factor(year) + factor(region) + factor(party_sub) | 0 | region, data=vrn_b, exactDOF=TRUE)

est3<-felm(natecon_avg ~ business_win+ natecon_year0 + totalexpend_log_year0 |  factor(unit_type)| 0 | region, data=subset(vrn_b, abs(businesswinmargin)<.03), exactDOF=TRUE)

est4<-with(vrn_b, rdrobust(y=natecon_avg , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      natecon_year0,
                                      totalexpend_log_year0
                           ), all=TRUE, kernel='uni',p=1,cluster=region,h=.05))

est5<-with(vrn_b, rdrobust(y=natecon_avg , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      natecon_year0,
                                      totalexpend_log_year0
                           ), all=TRUE, kernel='uni',p=1,cluster=region))

est6<-with(vrn_b, rdrobust(y=natecon_avg , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      factor(year),
                                      factor(region),
                                      factor(party_sub),
                                      natecon_year0,
                                      totalexpend_log_year0,
                                      transfers_level_year0,
                                      turnout,
                                      voterlist_log,
                                      numcandsperelection_log,
                                      incumbent), all=TRUE, kernel='uni',p=1,cluster=region))

####### Panel B

est1<-felm(education_avg ~ business_win+ education_year0 + totalexpend_log_year0 |  factor(unit_type) | 0 | region, data=vrn_b, exactDOF=TRUE)

est2<-felm(education_avg ~ business_win + education_year0 + totalexpend_log_year0  + transfers_level_year0 + turnout + voterlist_log + numcandsperelection_log + incumbent |  factor(unit_type) + factor(year) + factor(region) + factor(party_sub) | 0 | region, data=vrn_b, exactDOF=TRUE)

est3<-felm(education_avg ~ business_win+ education_year0 + totalexpend_log_year0 |  factor(unit_type)| 0 | region, data=subset(vrn_b, abs(businesswinmargin)<.03), exactDOF=TRUE)

est4<-with(vrn_b, rdrobust(y=education_avg , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      education_year0,
                                      totalexpend_log_year0
                           ), all=TRUE, kernel='uni',p=1,cluster=region,h=.05))

est5<-with(vrn_b, rdrobust(y=education_avg , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      education_year0,
                                      totalexpend_log_year0
                           ), all=TRUE, kernel='uni',p=1,cluster=region))

est6<-with(vrn_b, rdrobust(y=education_avg , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      factor(year),
                                      factor(region),
                                      factor(party_sub),
                                      education_year0,
                                      totalexpend_log_year0,
                                      transfers_level_year0,
                                      turnout,
                                      voterlist_log,
                                      numcandsperelection_log,
                                      incumbent), all=TRUE, kernel='uni',p=1,cluster=region))

#### Panel C

est1<-felm(health_avg ~ business_win+ health_year0 + totalexpend_log_year0 |  factor(unit_type) | 0 | region, data=vrn_b, exactDOF=TRUE)

est2<-felm(health_avg ~ business_win + health_year0 + totalexpend_log_year0  + transfers_level_year0 + turnout + voterlist_log + numcandsperelection_log + incumbent |  factor(unit_type) + factor(year) + factor(region) + factor(party_sub) | 0 | region, data=vrn_b, exactDOF=TRUE)

est3<-felm(health_avg ~ business_win+ health_year0 + totalexpend_log_year0 |  factor(unit_type)| 0 | region, data=subset(vrn_b, abs(businesswinmargin)<.03), exactDOF=TRUE)

est4<-with(vrn_b, rdrobust(y=health_avg , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      health_year0,
                                      totalexpend_log_year0
                           ), all=TRUE, kernel='uni',p=1,cluster=region,h=.05))

est5<-with(vrn_b, rdrobust(y=health_avg , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      health_year0,
                                      totalexpend_log_year0
                           ), all=TRUE, kernel='uni',p=1,cluster=region))

est6<-with(vrn_b, rdrobust(y=health_avg , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      factor(year),
                                      factor(region),
                                      factor(party_sub),
                                      health_year0,
                                      totalexpend_log_year0,
                                      transfers_level_year0,
                                      turnout,
                                      voterlist_log,
                                      numcandsperelection_log,
                                      incumbent), all=TRUE, kernel='uni',p=1,cluster=region))

#### Panel D

est1<-felm(totalexpend_avg_log ~ business_win+ totalexpend_log_year0 |  factor(unit_type) | 0 | region, data=vrn_b, exactDOF=TRUE)

est2<-felm(totalexpend_avg_log ~ business_win + totalexpend_log_year0  + transfers_level_year0 + turnout + voterlist_log + numcandsperelection_log + incumbent |  factor(unit_type) + factor(year) + factor(region) + factor(party_sub) | 0 | region, data=vrn_b, exactDOF=TRUE)

est3<-felm(totalexpend_avg_log ~ business_win+ totalexpend_log_year0 |  factor(unit_type)| 0 | region, data=subset(vrn_b, abs(businesswinmargin)<.03), exactDOF=TRUE)

est4<-with(vrn_b, rdrobust(y=totalexpend_avg_log , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      totalexpend_log_year0
                           ), all=TRUE, kernel='uni',p=1,cluster=region,h=.05))

est5<-with(vrn_b, rdrobust(y=totalexpend_avg_log , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      totalexpend_log_year0
                           ), all=TRUE, kernel='uni',p=1,cluster=region))

est6<-with(vrn_b, rdrobust(y=totalexpend_avg_log , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      factor(year),
                                      factor(region),
                                      factor(party_sub),
                                      totalexpend_log_year0,
                                      transfers_level_year0,
                                      turnout,
                                      voterlist_log,
                                      numcandsperelection_log,
                                      incumbent), all=TRUE, kernel='uni',p=1,cluster=region))

###########################################
########     TABLE   4      ################
###########################################

####### Panel A

est1<-felm(deficit_perc_avg ~ business_win+ deficit_perc_year0 + totalexpend_log_year0 |  factor(unit_type) | 0 | region, data=vrn_b, exactDOF=TRUE)

est2<-felm(deficit_perc_avg ~ business_win + deficit_perc_year0 + totalexpend_log_year0  + transfers_level_year0 + turnout + voterlist_log + numcandsperelection_log + incumbent |  factor(unit_type) + factor(year) + factor(region) + factor(party_sub) | 0 | region, data=vrn_b, exactDOF=TRUE)

est3<-felm(deficit_perc_avg ~ business_win+ deficit_perc_year0 + totalexpend_log_year0 |  factor(unit_type)| 0 | region, data=subset(vrn_b, abs(businesswinmargin)<.03), exactDOF=TRUE)

est4<-with(vrn_b, rdrobust(y=deficit_perc_avg , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      deficit_perc_year0,
                                      totalexpend_log_year0
                           ), all=TRUE, kernel='uni',p=1,cluster=region,h=.05))

est5<-with(vrn_b, rdrobust(y=deficit_perc_avg , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      deficit_perc_year0,
                                      totalexpend_log_year0
                           ), all=TRUE, kernel='uni',p=1,cluster=region))

est6<-with(vrn_b, rdrobust(y=deficit_perc_avg , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      factor(year),
                                      factor(region),
                                      factor(party_sub),
                                      deficit_perc_year0,
                                      totalexpend_log_year0,
                                      transfers_level_year0,
                                      turnout,
                                      voterlist_log,
                                      numcandsperelection_log,
                                      incumbent), all=TRUE, kernel='uni',p=1,cluster=region))

##### Panel B

est1<-felm(competitive_all ~ business_win + totalexpend_log_year0 |  factor(unit_type) | 0 | region, data=vrn_b, exactDOF=TRUE)

est2<-felm(competitive_all ~ business_win  + totalexpend_log_year0  + transfers_level_year0 + turnout + voterlist_log + numcandsperelection_log + incumbent |  factor(unit_type) + factor(year) + factor(region) + factor(party_sub) | 0 | region, data=vrn_b, exactDOF=TRUE)

est3<-felm(competitive_all ~ business_win + totalexpend_log_year0 |  factor(unit_type)| 0 | region, data=subset(vrn_b, abs(businesswinmargin)<.03), exactDOF=TRUE)

est4<-with(vrn_b, rdrobust(y=competitive_all , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      
                                      totalexpend_log_year0
                           ), all=TRUE, kernel='uni',p=1,cluster=region,h=.05))

est5<-with(vrn_b, rdrobust(y=competitive_all , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      
                                      totalexpend_log_year0
                           ), all=TRUE, kernel='uni',p=1,cluster=region))

est6<-with(vrn_b, rdrobust(y=competitive_all , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      factor(year),
                                      factor(region),
                                      factor(party_sub),
                                      
                                      totalexpend_log_year0,
                                      transfers_level_year0,
                                      turnout,
                                      voterlist_log,
                                      numcandsperelection_log,
                                      incumbent), all=TRUE, kernel='uni',p=1,cluster=region))

##### Panel C

est1<-felm(competitive_construction ~ business_win + totalexpend_log_year0 |  factor(unit_type) | 0 | region, data=vrn_b, exactDOF=TRUE)

est2<-felm(competitive_construction ~ business_win  + totalexpend_log_year0  + transfers_level_year0 + turnout + voterlist_log + numcandsperelection_log + incumbent |  factor(unit_type) + factor(year) + factor(region) + factor(party_sub) | 0 | region, data=vrn_b, exactDOF=TRUE)

est3<-felm(competitive_construction ~ business_win + totalexpend_log_year0 |  factor(unit_type)| 0 | region, data=subset(vrn_b, abs(businesswinmargin)<.03), exactDOF=TRUE)

est4<-with(vrn_b, rdrobust(y=competitive_construction , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      
                                      totalexpend_log_year0
                           ), all=TRUE, kernel='uni',p=1,cluster=region,h=.05))

est5<-with(vrn_b, rdrobust(y=competitive_construction , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      
                                      totalexpend_log_year0
                           ), all=TRUE, kernel='uni',p=1,cluster=region))

est6<-with(vrn_b, rdrobust(y=competitive_construction , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      factor(year),
                                      factor(region),
                                      factor(party_sub),
                                      
                                      totalexpend_log_year0,
                                      transfers_level_year0,
                                      turnout,
                                      voterlist_log,
                                      numcandsperelection_log,
                                      incumbent), all=TRUE, kernel='uni',p=1,cluster=region))

###########################################
########     FIGURE 1      ################
###########################################

plot.rdd <- function(varname,range,breaksize,data,label,title,newspan) {
  g=data
  g$variable <- g[,varname] 
  breaks <- c(sort(-seq(0, range, by = breaksize)[-1]), seq(0, range, by = breaksize))
  data.span  <- g[abs(g $businesswinmargin) <= max(breaks), ]
  data.span$bins <- cut(data.span$businesswinmargin, breaks, 
                        include.lowest = TRUE, right = FALSE)
  
  data.span.plot <- as.data.frame(cbind(tapply(data.span$variable, data.span$bins, mean, na.rm = TRUE)))
  data.span.N <- as.data.frame(cbind(tapply(data.span$variable, data.span$bins, length)))
  print(mean(data.span.N$V1))
  colnames(data.span.plot) <- c("variable")
  data.span.plot$runvar <- head(breaks, -1)
  data.span.plot$runvar = data.span.plot$runvar + breaksize/2
  p3 <- ggplot(data = subset(data.span.plot), mapping=aes(x = runvar, y = variable)) + geom_point() + stat_smooth(data = subset(data.span, businesswinmargin <=0), aes(x = businesswinmargin, y= variable),method = lm,span=newspan,level=.95) + stat_smooth(data = subset(data.span, businesswinmargin>0), aes(x = businesswinmargin, y= variable),method = lm,span=newspan) + ggtitle(title)+ xlab("Businessperson Mayor Margin of Victory")+ ylab(paste(label,"",sep=""))+ theme(axis.text = element_text(vjust=1,size=11),plot.title = element_text(hjust = 0.5))+theme_bw()+ theme(plot.title = element_text(hjust = 0.5),axis.title = element_text(size=15))
  p3
  return(p3) 
}

vrn_b_resid<-subset(vrn_b, is.na(deficit_perc_avg)==FALSE & is.na(deficit_perc_year0)==FALSE &is.na(totalexpend_log_year0)==FALSE)
est<-felm(deficit_perc_avg ~ deficit_perc_year0 + totalexpend_log_year0 |  factor(unit_type) | 0 | 0, data=vrn_b_resid, exactDOF=TRUE)
vrn_b_resid$plot_residuals<-est$residuals

p1<- plot.rdd("plot_residuals",range=.2,breaksize=.015,data=as.data.frame(vrn_b_resid), label="Deficit", title="(b)",newspan=.8)

vrn_b_resid<-subset(vrn_b, is.na(competitive_all)==FALSE & is.na(totalexpend_log_year0)==FALSE)
est<-felm(competitive_all ~ totalexpend_log_year0 |  factor(unit_type) | 0 | 0, data=vrn_b_resid, exactDOF=TRUE)
vrn_b_resid$plot_residuals<-est$residuals

p2<- plot.rdd("plot_residuals",range=.2,breaksize=.015,data=as.data.frame(vrn_b_resid), label="Competitive Procurement\nAll Purchases", title="(b)",newspan=.8)

vrn_b_resid<-subset(vrn_b, is.na(competitive_construction)==FALSE & is.na(totalexpend_log_year0)==FALSE)
est<-felm(competitive_construction ~ totalexpend_log_year0 |  factor(unit_type)| 0 | 0, data=vrn_b_resid, exactDOF=TRUE)
vrn_b_resid$plot_residuals<-est$residuals

p3<- plot.rdd("plot_residuals",range=.2,breaksize=.015,data=as.data.frame(vrn_b_resid), label="Competitive Procurement\nConstruction Purchases", title="(b)",newspan=.8)

vrn_b_resid<-subset(vrn_b, is.na(natecon_avg)==FALSE & is.na(natecon_year0)==FALSE &is.na(totalexpend_log_year0)==FALSE)
est<-felm(natecon_avg ~ natecon_year0 + totalexpend_log_year0 |  factor(unit_type) | 0 | 0, data=vrn_b_resid, exactDOF=TRUE)
vrn_b_resid$plot_residuals<-est$residuals

p4<- plot.rdd("plot_residuals",range=.2,breaksize=.015,data=as.data.frame(vrn_b_resid), label="Economic Infrastructure (%)", title="(d)",newspan=.8)

vrn_b_resid<-subset(vrn_b, is.na(education_avg)==FALSE & is.na(education_year0)==FALSE &is.na(totalexpend_log_year0)==FALSE)
est<-felm(education_avg ~ education_year0 + totalexpend_log_year0 |  factor(unit_type) | 0 | 0, data=vrn_b_resid, exactDOF=TRUE)
vrn_b_resid$plot_residuals<-est$residuals

p5<- plot.rdd("plot_residuals",range=.2,breaksize=.015,data=as.data.frame(vrn_b_resid), label="Education (%)", title="(e)",newspan=.8)

vrn_b_resid<-subset(vrn_b, is.na(health_avg)==FALSE & is.na(health_year0)==FALSE &is.na(totalexpend_log_year0)==FALSE)
est<-felm(health_avg ~ health_year0 + totalexpend_log_year0 |  factor(unit_type) | 0 | 0, data=vrn_b_resid, exactDOF=TRUE)
vrn_b_resid$plot_residuals<-est$residuals

p6<- plot.rdd("plot_residuals",range=.2,breaksize=.015,data=as.data.frame(vrn_b_resid), label="Health Care (%)", title="(f)",newspan=.8)

vrn_b_resid<-subset(vrn_b, is.na(totalexpend_avg_log)==FALSE &is.na(totalexpend_log_year0)==FALSE)
est<-felm(totalexpend_avg_log ~ totalexpend_log_year0 + totalexpend_log_year0 |  factor(unit_type) | 0 | 0, data=vrn_b_resid, exactDOF=TRUE)
vrn_b_resid$plot_residuals<-est$residuals

p7<- plot.rdd("plot_residuals",range=.2,breaksize=.015,data=as.data.frame(vrn_b_resid), label="Total Expenditures (log)", title="(a)",newspan=.8)








#########################################################
#############     SUPPLEMENTARY APPENDIX           ###################
#########################################################

###########################################
########     TABLE   A2      ################
###########################################

vrn_b$munirayon<-ifelse(vrn_b$unit_type=="муниципальный район",1,0)
vrn_b$gorpos<-ifelse(vrn_b$unit_type=="городское поселение",1,0)
vrn_b$gorokr<-ifelse(vrn_b$unit_type=="городской округ",1,0)
vrn_b$selpos<-ifelse(vrn_b$unit_type=="сельское поселение",1,0)

vrn_b<-as.data.table(vrn_b)
vrn_summary<-vrn_b[,list(business_win,
                         businesswinmargin,
                         munirayon,
                         gorpos,
                         selpos,
                         gorokr,
                         averageterritory,
                         voterlist_log,
                         turnout,
                         numcandsperelection,
                         election_urwin,
                         incumbent,
                         log_age,
                         female,
                         deficit_perc_avg,
                         competitive_all,
                         competitive_construction,
                         natecon_avg,
                         education_avg,
                         health_avg,
                         totalexpend_avg_log,
                         commongov_avg,
                         transfers_avg,
                         totalexpend_log_year0,
                         transfers_level_year0,
                         education_year0,
                         health_year0,
                         natecon_year0,
                         commongov_year0,
                         deficit_perc_year0)]


vrntable<-mod_stargazer(stargazer(vrn_summary,covariate.labels=c("Businessperson Mayor","Businessperson Vote Margin","Municipal Rayon","City Settlement","Rural Settlement","City District","Geographic Size (log)","Number of Voters on List (log)","Turnout","Number of Candidates","United Russia Mayor","Businessperson Incumbent","Businessperson Age (log)","Businessperson Female",
                                                                 "Deficit (Expenditures/Revenue)",
                                                                 "Competitive Procurement - All Purchases",
                                                                 "Competitive Procurement - Construction Purchases", 
                                                                 "Economic Infrastructure (as \\% of expenditures)",
                                                                 "Education (as \\% of expenditures)","Health (as \\% of expenditures)",
                                                                 "Total Expenditures (log)","Government Administration (as \\% of expenditures)",
                                                                 "Transfers (as \\% of revenue)","Total Expenditures (log)","Transfers (as \\% of revenue)","Education (as \\% of expenditures)","Health (as \\% of expenditures)","Economic Infrastructure (as \\% of expenditures)","Government Administration (as \\% of expenditures)","Deficit (Expenditures/Revenue)"),header=FALSE,digits=3))


###########################################
########     FIGURE   B1      ################
###########################################

exp1<-ggplot(data=subset(vrn_b, munirayon==1), aes(totalexpend_avg_log))+geom_histogram(col="royalblue4",fill="royalblue4",alpha=.2, binwidth = .25) + labs(x="Total Expenditures (log)", y="Count")+xlim(6,18)+ theme(
  plot.title = element_text(size=12, face="plain"),
  axis.title = element_text(size=10),axis.text = element_text(size=10))

exp2<-ggplot(data=subset(vrn_b, gorokr==1), aes(totalexpend_avg_log))+geom_histogram(col="royalblue4",fill="royalblue4",alpha=.2, binwidth = .25) + labs(x="Total Expenditures (log)", y="Count")+xlim(6,18)+ theme(
  plot.title = element_text(size=12, face="plain"),
  axis.title = element_text(size=10),axis.text = element_text(size=10))

exp3<-ggplot(data=subset(vrn_b, gorpos==1), aes(totalexpend_avg_log))+geom_histogram(col="royalblue4",fill="royalblue4",alpha=.2, binwidth = .25) + labs(x="Total Expenditures (log)", y="Count")+xlim(6,18)+ theme(
  plot.title = element_text(size=12, face="plain"),
  axis.title = element_text(size=10),axis.text = element_text(size=10))

exp4<-ggplot(data=subset(vrn_b, selpos==1), aes(totalexpend_avg_log))+geom_histogram(col="royalblue4",fill="royalblue4",alpha=.2, binwidth = .25) + labs(x="Total Expenditures (log)", y="Count")+xlim(6,18)+ theme(
  plot.title = element_text(size=12, face="plain"),
  axis.title = element_text(size=10),axis.text = element_text(size=10))

def1<-ggplot(data=subset(vrn_b, munirayon==1), aes(deficit_perc_avg))+geom_histogram(col="royalblue4",fill="royalblue4",alpha=.2, binwidth = .025) + labs(x="Budget Deficit (%)", y="Count")+xlim(.5,1.5)+ theme(
  plot.title = element_text(size=12, face="plain"),
  axis.title = element_text(size=10),axis.text = element_text(size=10))

def2<-ggplot(data=subset(vrn_b, gorokr==1), aes(deficit_perc_avg))+geom_histogram(col="royalblue4",fill="royalblue4",alpha=.2, binwidth = .025) + labs(x="Budget Deficit (%)", y="Count")+xlim(.5,1.5)+ theme(
  plot.title = element_text(size=12, face="plain"),
  axis.title = element_text(size=10),axis.text = element_text(size=10))

def3<-ggplot(data=subset(vrn_b, gorpos==1), aes(deficit_perc_avg))+geom_histogram(col="royalblue4",fill="royalblue4",alpha=.2, binwidth = .025) + labs(x="Budget Deficit (%)", y="Count")+xlim(.5,1.5)+ theme(
  plot.title = element_text(size=12, face="plain"),
  axis.title = element_text(size=10),axis.text = element_text(size=10))

def4<-ggplot(data=subset(vrn_b, selpos==1), aes(deficit_perc_avg))+geom_histogram(col="royalblue4",fill="royalblue4",alpha=.2, binwidth = .025) + labs(x="Budget Deficit (%)", y="Count")+xlim(.5,1.5)+ theme(
  plot.title = element_text(size=12, face="plain"),
  axis.title = element_text(size=10),axis.text = element_text(size=10))

comp1<-ggplot(data=subset(vrn_b, munirayon==1), aes(competitive_all))+geom_histogram(col="royalblue4",fill="royalblue4",alpha=.2, binwidth = .025) + labs(x="Competitive Procurement (all)", y="Count")+xlim(0,1)+ theme(
  plot.title = element_text(size=12, face="plain"),
  axis.title = element_text(size=10),axis.text = element_text(size=10))

comp2<-ggplot(data=subset(vrn_b, gorokr==1), aes(competitive_all))+geom_histogram(col="royalblue4",fill="royalblue4",alpha=.2, binwidth = .025) + labs(x="Competitive Procurement (all)", y="Count")+xlim(0,1)+ theme(
  plot.title = element_text(size=12, face="plain"),
  axis.title = element_text(size=10),axis.text = element_text(size=10))

comp3<-ggplot(data=subset(vrn_b, gorpos==1), aes(competitive_all))+geom_histogram(col="royalblue4",fill="royalblue4",alpha=.2, binwidth = .025) + labs(x="Competitive Procurement (all)", y="Count")+xlim(0,1)+ theme(
  plot.title = element_text(size=12, face="plain"),
  axis.title = element_text(size=10),axis.text = element_text(size=10))

comp4<-ggplot(data=subset(vrn_b, selpos==1), aes(competitive_all))+geom_histogram(col="royalblue4",fill="royalblue4",alpha=.2, binwidth = .025) + labs(x="Competitive Procurement (all)", y="Count")+xlim(0,1)+ylim(0,100)+ theme(
  plot.title = element_text(size=12, face="plain"),
  axis.title = element_text(size=10),axis.text = element_text(size=10))


def1<-ggplot(data=subset(vrn_b, munirayon==1), aes(natecon_avg))+geom_histogram(col="royalblue4",fill="royalblue4",alpha=.2, binwidth = .025) + labs(x="Economic Expenditures (%)", y="Count")+xlim(0,1)+ theme(
  plot.title = element_text(size=12, face="plain"),
  axis.title = element_text(size=10),axis.text = element_text(size=10))

def2<-ggplot(data=subset(vrn_b, gorokr==1), aes(natecon_avg))+geom_histogram(col="royalblue4",fill="royalblue4",alpha=.2, binwidth = .025) + labs(x="Economic Expenditures (%)", y="Count")+ylim(0,50)+xlim(0,1)+ theme(
  plot.title = element_text(size=12, face="plain"),
  axis.title = element_text(size=10),axis.text = element_text(size=10))

def3<-ggplot(data=subset(vrn_b, gorpos==1), aes(natecon_avg))+geom_histogram(col="royalblue4",fill="royalblue4",alpha=.2, binwidth = .025) + labs(x="Economic Expenditures (%)", y="Count")+xlim(0,1)+ylim(0,50)+ theme(
  plot.title = element_text(size=12, face="plain"),
  axis.title = element_text(size=10),axis.text = element_text(size=10))

def4<-ggplot(data=subset(vrn_b, selpos==1), aes(natecon_avg))+geom_histogram(col="royalblue4",fill="royalblue4",alpha=.2, binwidth = .025) + labs(x="Economic Expenditures (%)", y="Count")+xlim(0,1)+ theme(
  plot.title = element_text(size=12, face="plain"),
  axis.title = element_text(size=10),axis.text = element_text(size=10))



exp1<-ggplot(data=subset(vrn_b, munirayon==1), aes(education_avg))+geom_histogram(col="royalblue4",fill="royalblue4",alpha=.2, binwidth = .025) + labs(x="Education Expenditures (%)", y="Count")+xlim(0,1)+ylim(0,50)+ theme(
  plot.title = element_text(size=12, face="plain"),
  axis.title = element_text(size=10),axis.text = element_text(size=10))

exp2<-ggplot(data=subset(vrn_b, gorokr==1), aes(education_avg))+geom_histogram(col="royalblue4",fill="royalblue4",alpha=.2, binwidth = .025) + labs(x="Education Expenditures (%)", y="Count")+xlim(0,1)+ylim(0,50)+ theme(
  plot.title = element_text(size=12, face="plain"),
  axis.title = element_text(size=10),axis.text = element_text(size=10))

exp3<-ggplot(data=subset(vrn_b, gorpos==1), aes(education_avg))+geom_histogram(col="royalblue4",fill="royalblue4",alpha=.2, binwidth = .025) + labs(x="Education Expenditures (%)", y="Count")+xlim(0,1)+ theme(
  plot.title = element_text(size=12, face="plain"),
  axis.title = element_text(size=10),axis.text = element_text(size=10))+ylim(0,100)

exp4<-ggplot(data=subset(vrn_b, selpos==1), aes(education_avg))+geom_histogram(col="royalblue4",fill="royalblue4",alpha=.2, binwidth = .025) + labs(x="Education Expenditures (%)", y="Count")+xlim(0,1)+ theme(
  plot.title = element_text(size=12, face="plain"),
  axis.title = element_text(size=10),axis.text = element_text(size=10))+ylim(0,100)



health1<-ggplot(data=subset(vrn_b, munirayon==1), aes(health_avg))+geom_histogram(col="royalblue4",fill="royalblue4",alpha=.2, binwidth = .025) + labs(x="Health Expenditures (%)", y="Count")+xlim(0,1)+ theme(
  plot.title = element_text(size=12, face="plain"),
  axis.title = element_text(size=10),axis.text = element_text(size=10))

health2<-ggplot(data=subset(vrn_b, gorokr==1), aes(health_avg))+geom_histogram(col="royalblue4",fill="royalblue4",alpha=.2, binwidth = .025) + labs(x="Health Expenditures (%)", y="Count")+xlim(0,1)+ylim(0,50)+ theme(
  plot.title = element_text(size=12, face="plain"),
  axis.title = element_text(size=10),axis.text = element_text(size=10))

health3<-ggplot(data=subset(vrn_b, gorpos==1), aes(health_avg))+geom_histogram(col="royalblue4",fill="royalblue4",alpha=.2, binwidth = .025) + labs(x="Health Expenditures (%)", y="Count")+xlim(0,1)+ theme(
  plot.title = element_text(size=12, face="plain"),
  axis.title = element_text(size=10),axis.text = element_text(size=10))+ylim(0,100)

health4<-ggplot(data=subset(vrn_b, selpos==1), aes(health_avg))+geom_histogram(col="royalblue4",fill="royalblue4",alpha=.2, binwidth = .025) + labs(x="Health Expenditures (%)", y="Count")+xlim(0,1)+ theme(
  plot.title = element_text(size=12, face="plain"),
  axis.title = element_text(size=10),axis.text = element_text(size=10))+ylim(0,250)



###########################################
########     FIGURE   C1      ################
###########################################
binsize = .01

se <- DCdensity(vrn_b_year0$businesswinmargin,0, bin= binsize, ext.out=TRUE)$se
t <- DCdensity(vrn_b_year0$businesswinmargin,0, bin= binsize, ext.out=TRUE)$theta
p <- DCdensity(vrn_b_year0$businesswinmargin,0, bin= binsize, ext.out=TRUE)$p
t <- round(t, digits = 3)
se <- round(se, digits = 3)
p <- round(p, digits = 3)
abline(v=0,lty=3, col = 'red')
mtext(paste("Difference=",t), side=1, line=2, cex=.8)
mtext(paste("Standard Error=",se), side=1, line=3, cex=.8)
mtext(paste("p-value=",p), side=1, line=4, cex=.8)


###########################################
########     TABLE C1 / FIGURE   C2      ################
###########################################

#### These loops save both the t-statistics and the regression objects for Table C1 and Figure C2

### List of placebo variables to run regressions on
dvs<-c("totalexpend_log_year0", 
       "education_year0",
       "zhkx_year0", 
       "culture_year0",
       "natecon_year0",
       "commongov_year0", 
       "health_year0",
       "security_year0",
       "totalincome_log_year0",
       "transfers_level_year0", 
       "deficit_perc_year0",
       "perc_ch_income_year0",
       "housing_new_log", 
       "voterlist_log",
       "numcandsperelection_log",
       "turnout",
       "ur",
       "incumbent",
       "log_age",
       "female")

numv<-length(dvs)

cm_3_t = rep(0,length(dvs))
cm_5_t = rep(0,length(dvs))
ll_05_t = rep(0,length(dvs))
ll_optb_t = rep(0,length(dvs))
opts1 = rep(0,length(dvs))

ll_5_m_c1 = list()
ll_10_m_c1 = list()
ll_10_n_c1 = list()

vrn_b_year0 <-as.data.frame(vrn_b_year0)

for(i in 1:length(dvs)){
  print(i)
  cm_3_t[i] = felm(formula(paste(substitute(i,list(i = as.name(dvs[i]))), "~ business_win |  factor(unit_type) | 0 | region",sep="")), data=subset(vrn_b_year0, abs(businesswinmargin)<.03), exactDOF=TRUE)$ctval["business_win"]
  
  cm_5_t[i] = felm(formula(paste(substitute(i,list(i = as.name(dvs[i]))), "~ business_win |  factor(unit_type) | 0 | region",sep="")), data=subset(vrn_b_year0, abs(businesswinmargin)<.05), exactDOF=TRUE)$ctval["business_win"]
  
  ll_05_t[i] = abs(with(vrn_b_year0, rdrobust(y=get(paste(substitute(i,list(i = as.name(dvs[i]))))) , x=businesswinmargin, bwselect="mserd",
                                                          covs=cbind(factor(unit_type)), all=TRUE, kernel='uni',p=1,cluster=region,h=.1))$z[3])
  
  ll_optb_t[i] = abs(with(vrn_b_year0, rdrobust(y=get(paste(substitute(i,list(i = as.name(dvs[i]))))) , x=businesswinmargin, bwselect="mserd",
                                                            covs=cbind(factor(unit_type)), all=TRUE, kernel='uni',p=1,cluster=region))$z[3])
  
  ## Extracting Coefs for table - 
  ll_5_m_c1[[i]] = felm(formula(paste(substitute(i,list(i = as.name(dvs[i]))), "~ business_win |  factor(unit_type) | 0 | region",sep="")), data=subset(vrn_b_year0, abs(businesswinmargin)<.05), exactDOF=TRUE)
  
  ## Extracting Coefs for table - 
  ll_10_m_c1[[i]] = with(vrn_b_year0, rdrobust(y=get(paste(substitute(i,list(i = as.name(dvs[i]))))) , x=businesswinmargin, bwselect="mserd",
                                                           covs=cbind(factor(unit_type)), all=TRUE, kernel='uni',p=1,cluster=region))
  
  ll_10_n_c1[[i]] = felm(formula(paste(substitute(i,list(i = as.name(dvs[i]))), "~ business_win |  factor(unit_type) | 0 | region",sep="")), data=subset(vrn_b_year0, abs(businesswinmargin)<ll_10_m_c1[[i]]$bws[3]), exactDOF=TRUE)
  
}

set_t <- cbind(cm_3_t, cm_5_t, ll_05_t,ll_optb_t)


covs_labels <- c( "Total Expenditures (log)",
                  "Education (% of exp.)",
                  "Housing (% of exp.)",
                  "Culture (% of exp.)",
                  "Economic Infr. (% of exp.)",
                  "Government Oper. (% of exp.)",
                  "Health (% of exp.)",
                  "Security (% of exp.)",
                  "Total Revenue (log)",
                  "Transfers (as % of revenue)",
                  "Deficit (Expenditures/Revenue)",
                  "Change in Revenue (%)", 
                  "New Housing Construction (log)",
                  "Number of Voters on List (log)",
                  "Number of Candidates (log)",
                  "Turnout",
                  "Member of United Russia",
                  "Incumbent",
                  "Age",
                  "Female")

tstat<-as.data.frame(set_t)
tstat$Variable <- covs_labels
mdata <- melt(tstat, id=c("Variable"))
tstat1<-head(mdata, -1*(numv))

tstat1<-mdata
levels(tstat1$variable) <- c("Bandwidth = 3%",'Bandwidth = 5%','Bandwidth = 10%',"Bandwidth = Optimal")
tstat1$value<-abs(as.numeric(as.character(tstat1$value)))

Estimator = c(rep("Close Margin", numv*2),rep("Local Linear", numv*2))
tstat1$Estimator <- c(rep("Close Margin", numv*2),rep("Local Linear", numv*2))
names(tstat1)[2] <- "Sample"
ggplot(tstat1, aes(x = value, y = Variable, shape = Sample, color = Sample, size=Sample)) + geom_point() + geom_vline(xintercept = 2, lty = 2, size = .2) + facet_grid(.~Estimator) + theme_bw() + scale_x_continuous("T Statistic") + scale_color_brewer(palette = "Set1") + theme(axis.text = element_text(size=12),legend.text=element_text(size=12),legend.title=element_text(size=12),strip.text.x = element_text(size = 12),axis.title = element_text(size=12))+scale_size_manual(values=c(3,3,3,3))




###########################################
########     TABLE   D1      ################
###########################################

### Just take municipalities with elections
vrn_d<-subset(exp_full, is.na(numbuscands)==FALSE)

vrn_d$hasbusinesscandidate<-ifelse(vrn_d$numbuscands>0,1,0)
vrn_d$unit_type<-as.character(vrn_d$unit_type)
vrn_d<-subset(vrn_d, unit_type!="")

vrn_d$education<-vrn_d$eduspend/vrn_d$totalexpend
vrn_d$zhkx<-vrn_d$zhkxspend/vrn_d$totalexpend
vrn_d$culture<-vrn_d$culturespend/vrn_d$totalexpend
vrn_d$natecon<-vrn_d$nateconspend/vrn_d$totalexpend
vrn_d$commongov<-vrn_d$commongovexpend/vrn_d$totalexpend
vrn_d$social<-vrn_d$socialspend/vrn_d$totalexpend
vrn_d$security<-vrn_d$securityspend/vrn_d$totalexpend
vrn_d$health<-vrn_d$totalhealthspend/vrn_d$totalexpend
vrn_d$housing_new_log<-log(vrn_d$housing_new+1)
vrn_d$transfers_level<-vrn_d$transfers/vrn_d$totalincome


vrn_d$economicdata_expend<-ifelse(is.na(vrn_d$totalexpend)==FALSE,1,0)
vrn_d$economicdata_natecon<-ifelse(is.na(vrn_d$natecon)==FALSE,1,0)
vrn_d$economicdata_commongov<-ifelse(is.na(vrn_d$commongov)==FALSE,1,0)
vrn_d$economicdata_education<-ifelse(is.na(vrn_d$education)==FALSE,1,0)
vrn_d$economicdata_zhkx<-ifelse(is.na(vrn_d$zhkx)==FALSE,1,0)
vrn_d$economicdata_health<-ifelse(is.na(vrn_d$health)==FALSE,1,0)
vrn_d$economicdata_culture<-ifelse(is.na(vrn_d$culture)==FALSE,1,0)
vrn_d$economicdata_contracts<-ifelse(is.na(vrn_d$total_contracts)==FALSE,1,0)

vrn_d$turnout<-vrn_d$electiontotalvotes/vrn_d$voterlist
vrn_d$electionyear<-vrn_d$year

est1 <- felm(hasbusinesscandidate ~ unit_type + log(voterlist) + ur + incumbent + factor(electionyear)| region | 0 | region, data=vrn_d)

est2 <- felm(hasbusinesscandidate ~ unit_type + log(voterlist) + ur  + incumbent +factor(electionyear) + totalexpend_log + natecon + commongov + transfers_level| region | 0 | region, data=vrn_d)

est3 <- felm(numbuscands ~ incumbent + unit_type + log(voterlist) +incumbent + factor(electionyear)| region | 0 | region, data=vrn_d)

est4 <- felm(numbuscands ~ unit_type + log(voterlist) + ur  + incumbent +factor(electionyear) + totalexpend_log + natecon + commongov + transfers_level| region | 0 | region, data=vrn_d)

est5 <- felm(businesswin ~ unit_type + log(voterlist) + ur  +incumbent + factor(electionyear)| region | 0 | region, data=vrn_d)

est6 <- felm(businesswin ~ unit_type + log(voterlist) + ur  + incumbent + factor(electionyear) + totalexpend_log + natecon + commongov + transfers_level | region | 0 | region, data=vrn_d)


###########################################
########     TABLE   D2      ################
###########################################


est1 <- felm(economicdata_expend ~ unit_type + log(voterlist) + ur + incumbent + hasbusinesscandidate + turnout + factor(electionyear)| region | 0 | region, data=vrn_d)

est2 <- felm(economicdata_natecon ~ unit_type + log(voterlist) + ur + incumbent + hasbusinesscandidate + turnout + factor(electionyear)| region | 0 | region, data=vrn_d)

est3 <- felm(economicdata_commongov ~ unit_type + log(voterlist) + ur + incumbent + hasbusinesscandidate + turnout + factor(electionyear)| region | 0 | region, data=vrn_d)

est4 <- felm(economicdata_education ~ unit_type + log(voterlist) + ur + incumbent + hasbusinesscandidate + turnout + factor(electionyear)| region | 0 | region, data=vrn_d)

est5 <- felm(economicdata_zhkx ~ unit_type + log(voterlist) + ur + incumbent + hasbusinesscandidate+ turnout  + factor(electionyear)| region | 0 | region, data=vrn_d)

est6 <- felm(economicdata_culture ~ unit_type + log(voterlist) + ur + incumbent + hasbusinesscandidate + turnout + factor(electionyear)| region | 0 | region, data=vrn_d)

est7 <- felm(economicdata_health ~ unit_type + log(voterlist) + ur + incumbent + hasbusinesscandidate+ turnout  + factor(electionyear)| region | 0 | region, data=vrn_d)

est8 <- felm(economicdata_contracts ~ unit_type + log(voterlist) + ur + incumbent + hasbusinesscandidate+ turnout  + factor(electionyear)| region | 0 | region, data=subset(vrn_d, electionyear>2010))

###########################################
########     TABLE   D3      ################
###########################################

### Average out over years
mundata_mun<-exp_full[,list(totalexpend=mean(totalexpend,na.rm=TRUE),
                            totalincome=mean(totalincome,na.rm=TRUE),
                            transfers=mean(transfers,na.rm=TRUE),
                            
                            allpopulation=mean(allpopulation,na.rm=TRUE),
                            totalterritory=mean(totalterritory,na.rm=TRUE),
                            numyears=uniqueN(year),
                            numelections=sum(haselection)), by=c("oktmo","unit_type")
                      ]
mundata_mun$haselection<-ifelse(mundata_mun$numelections>0,1,0)

mundata_mun<-subset(mundata_mun, unit_type!="внутригородская территория города федерального значения")

mundata_mun$region<-str_sub(mundata_mun$oktmo,0,2)

mundata_mun$transfers_level<-mundata_mun$transfers/mundata_mun$totalincome
mundata_mun$deficit<-mundata_mun$totalexpend/mundata_mun$totalincome
mundata_mun$population_log<-log(mundata_mun$allpopulation)
mundata_mun$totalterritory_log<-log(mundata_mun$totalterritory)
mundata_mun$totalexpend_log<-log(mundata_mun$totalexpend)

mundata_mun[mundata_mun=="-Inf"]<-NA

est1 <- felm(haselection ~ unit_type + totalterritory_log +totalexpend_log + transfers_level + deficit | 0 | 0 | region, data=mundata_mun)

est2 <- felm(haselection ~ unit_type + totalterritory_log +totalexpend_log + transfers_level + deficit | factor(region) | 0 | region, data=mundata_mun)

est3 <- felm(haselection ~ unit_type + totalterritory_log +totalexpend_log + transfers_level + deficit + population_log | factor(region) | 0 | region, data=mundata_mun)




###########################################
########     TABLE   E1      ################
###########################################

####### Panel A

est1<-felm(large_deficit ~ business_win+ deficit_perc_year0 + totalexpend_log_year0 |  factor(unit_type) | 0 | region, data=vrn_b, exactDOF=TRUE)

est2<-felm(large_deficit ~ business_win + deficit_perc_year0 + totalexpend_log_year0  + transfers_level_year0 + turnout + voterlist_log + numcandsperelection_log + incumbent |  factor(unit_type) + factor(year) + factor(region) + factor(party_sub) | 0 | region, data=vrn_b, exactDOF=TRUE)

est3<-felm(large_deficit ~ business_win+ deficit_perc_year0 + totalexpend_log_year0 |  factor(unit_type)| 0 | region, data=subset(vrn_b, abs(businesswinmargin)<.03), exactDOF=TRUE)

est4<-with(vrn_b, rdrobust(y=large_deficit , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      deficit_perc_year0,
                                      totalexpend_log_year0
                           ), all=TRUE, kernel='uni',p=1,cluster=region,h=.05))

est5<-with(vrn_b, rdrobust(y=large_deficit , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      deficit_perc_year0,
                                      totalexpend_log_year0
                           ), all=TRUE, kernel='uni',p=1,cluster=region))

est6<-with(vrn_b, rdrobust(y=large_deficit , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      factor(year),
                                      factor(region),
                                      factor(party_sub),
                                      deficit_perc_year0,
                                      totalexpend_log_year0,
                                      transfers_level_year0,
                                      turnout,
                                      voterlist_log,
                                      numcandsperelection_log,
                                      incumbent), all=TRUE, kernel='uni',p=1,cluster=region))


#### Panel B

est1<-felm(transfers_avg ~ business_win+ transfers_level_year0 + totalexpend_log_year0 |  factor(unit_type) | 0 | region, data=vrn_b, exactDOF=TRUE)

est2<-felm(transfers_avg ~ business_win + transfers_level_year0 + totalexpend_log_year0  + transfers_level_year0 + turnout + voterlist_log + numcandsperelection_log + incumbent |  factor(unit_type) + factor(year) + factor(region) + factor(party_sub) | 0 | region, data=vrn_b, exactDOF=TRUE)

est3<-felm(transfers_avg ~ business_win+ transfers_level_year0 + totalexpend_log_year0 |  factor(unit_type)| 0 | region, data=subset(vrn_b, abs(businesswinmargin)<.03), exactDOF=TRUE)

est4<-with(vrn_b, rdrobust(y=transfers_avg , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      transfers_level_year0,
                                      totalexpend_log_year0
                           ), all=TRUE, kernel='uni',p=1,cluster=region,h=.05))

est5<-with(vrn_b, rdrobust(y=transfers_avg , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      transfers_level_year0,
                                      totalexpend_log_year0
                           ), all=TRUE, kernel='uni',p=1,cluster=region))

est6<-with(vrn_b, rdrobust(y=transfers_avg , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      factor(year),
                                      factor(region),
                                      factor(party_sub),
                                      transfers_level_year0,
                                      totalexpend_log_year0,
                                      turnout,
                                      voterlist_log,
                                      numcandsperelection_log,
                                      incumbent), all=TRUE, kernel='uni',p=1,cluster=region))


###########################################
########     TABLE   E2      ################
###########################################

est1<-with(vrn_b, rdrobust(y=natecon_avg , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      natecon_year0,
                                      totalexpend_log_year0
                           ), all=TRUE, kernel='uni',p=1,cluster=region))

est2<-with(vrn_b, rdrobust(y=education_avg , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      education_year0,
                                      totalexpend_log_year0
                           ), all=TRUE, kernel='uni',p=1,cluster=region))

est3<-with(vrn_b, rdrobust(y=health_avg , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      health_year0,
                                      totalexpend_log_year0
                           ), all=TRUE, kernel='uni',p=1,cluster=region))

est4<-with(vrn_b, rdrobust(y=zhkx_avg , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      zhkx_year0,
                                      totalexpend_log_year0
                           ), all=TRUE, kernel='uni',p=1,cluster=region))

est5<-with(vrn_b, rdrobust(y=social_avg , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      factor(region),
                                      factor(year),
                                      social_year0,
                                      totalexpend_log_year0
                           ), all=TRUE, kernel='uni',p=1,cluster=region))

est6<-with(vrn_b, rdrobust(y=commongov_avg , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      commongov_year0,
                                      totalexpend_log_year0
                           ), all=TRUE, kernel='uni',p=1,cluster=region))

est7<-with(vrn_b, rdrobust(y=culture_avg , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      culture_year0,
                                      totalexpend_log_year0
                           ), all=TRUE, kernel='uni',p=1,cluster=region))


###########################################
########     TABLE   E3      ################
###########################################

est1<-with(vrn_b, rdrobust(y=competitive_construction , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      
                                      totalexpend_log_year0
                           ), all=TRUE, kernel='uni',p=1,cluster=region))

est2<-with(vrn_b, rdrobust(y=competitive_construction , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      factor(year),
                                      factor(region),
                                      factor(party_sub),
                                      
                                      totalexpend_log_year0,
                                      transfers_level_year0,
                                      turnout,
                                      voterlist_log,
                                      numcandsperelection_log,
                                      incumbent), all=TRUE, kernel='uni',p=1,cluster=region))

est3<-with(vrn_b, rdrobust(y=competitive_cars , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      
                                      totalexpend_log_year0
                           ), all=TRUE, kernel='uni',p=1,cluster=region))

est4<-with(vrn_b, rdrobust(y=competitive_cars , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      factor(year),
                                      factor(region),
                                      factor(party_sub),
                                      
                                      totalexpend_log_year0,
                                      transfers_level_year0,
                                      turnout,
                                      voterlist_log,
                                      numcandsperelection_log,
                                      incumbent), all=TRUE, kernel='uni',p=1,cluster=region))

est5<-with(vrn_b, rdrobust(y=competitive_food , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      
                                      totalexpend_log_year0
                           ), all=TRUE, kernel='uni',p=1,cluster=region))

est6<-with(vrn_b, rdrobust(y=competitive_food , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      factor(year),
                                      factor(region),
                                      factor(party_sub),
                                      
                                      totalexpend_log_year0,
                                      transfers_level_year0,
                                      turnout,
                                      voterlist_log,
                                      numcandsperelection_log
                           ), all=TRUE, kernel='uni',p=1,cluster=region))

est7<-with(vrn_b, rdrobust(y=competitive_office , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      
                                      totalexpend_log_year0
                           ), all=TRUE, kernel='uni',p=1,cluster=region))

est8<-with(vrn_b, rdrobust(y=competitive_office , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      factor(year),
                                      factor(region),
                                      factor(party_sub),
                                      
                                      totalexpend_log_year0,
                                      transfers_level_year0,
                                      turnout,
                                      voterlist_log,
                                      numcandsperelection_log,
                                      incumbent), all=TRUE, kernel='uni',p=1,cluster=region))

est9<-with(vrn_b, rdrobust(y=competitive_furniture , x=businesswinmargin, bwselect="mserd",
                           covs=cbind(factor(unit_type),
                                      
                                      totalexpend_log_year0
                           ), all=TRUE, kernel='uni',p=1,cluster=region))

est10<-with(vrn_b, rdrobust(y=competitive_furniture , x=businesswinmargin, bwselect="mserd",
                            covs=cbind(factor(unit_type),
                                       factor(year),
                                       factor(region),
                                       factor(party_sub),
                                       
                                       totalexpend_log_year0,
                                       transfers_level_year0,
                                       turnout,
                                       voterlist_log,
                                       numcandsperelection_log,
                                       incumbent), all=TRUE, kernel='uni',p=1,cluster=region))



###########################################
########     TABLE   E4      ################
###########################################

#### Panel A

est1<-with(subset(vrn_b, reg_democracy<=30), rdrobust(y=natecon_avg , x=businesswinmargin, bwselect="mserd",
                                                      covs=cbind(factor(unit_type),
                                                                 natecon_year0,
                                                                 totalexpend_log_year0
                                                      ), all=TRUE, kernel='uni',p=1,cluster=region))

est2<-with(subset(vrn_b, reg_democracy>30), rdrobust(y=natecon_avg , x=businesswinmargin, bwselect="mserd",
                                                     covs=cbind(factor(unit_type),
                                                                natecon_year0,
                                                                totalexpend_log_year0
                                                     ), all=TRUE, kernel='uni',p=1,cluster=region))

est3<-with(subset(vrn_b, reg_democracy<=30), rdrobust(y=education_avg , x=businesswinmargin, bwselect="mserd",
                                                      covs=cbind(factor(unit_type),
                                                                 education_year0,
                                                                 totalexpend_log_year0
                                                      ), all=TRUE, kernel='uni',p=1,cluster=region))

est4<-with(subset(vrn_b, reg_democracy>30), rdrobust(y=education_avg , x=businesswinmargin, bwselect="mserd",
                                                     covs=cbind(factor(unit_type),
                                                                education_year0,
                                                                totalexpend_log_year0
                                                     ), all=TRUE, kernel='uni',p=1,cluster=region))

est5<-with(subset(vrn_b, reg_democracy<=30), rdrobust(y=health_avg , x=businesswinmargin, bwselect="mserd",
                                                      covs=cbind(factor(unit_type),
                                                                 health_year0,
                                                                 totalexpend_log_year0
                                                      ), all=TRUE, kernel='uni',p=1,cluster=region))

est6<-with(subset(vrn_b, reg_democracy>30), rdrobust(y=health_avg , x=businesswinmargin, bwselect="mserd",
                                                     covs=cbind(factor(unit_type),
                                                                health_year0,
                                                                totalexpend_log_year0
                                                     ), all=TRUE, kernel='uni',p=1,cluster=region))

est7<-with(subset(vrn_b, reg_democracy<=30), rdrobust(y=totalexpend_avg_log , x=businesswinmargin, bwselect="mserd",
                                                      covs=cbind(factor(unit_type),
                                                                 totalexpend_log_year0
                                                      ), all=TRUE, kernel='uni',p=1,cluster=region))

est8<-with(subset(vrn_b, reg_democracy>30), rdrobust(y=totalexpend_avg_log , x=businesswinmargin, bwselect="mserd",
                                                     covs=cbind(factor(unit_type),
                                                                totalexpend_log_year0
                                                     ), all=TRUE, kernel='uni',p=1,cluster=region))

### Panel B

est1<-with(subset(vrn_b, reg_democracy<=30), rdrobust(y=deficit_perc_avg , x=businesswinmargin, bwselect="mserd",
                                                      covs=cbind(factor(unit_type),
                                                                 deficit_perc_year0,
                                                                 totalexpend_log_year0
                                                      ), all=TRUE, kernel='uni',p=1,cluster=region))

est2<-with(subset(vrn_b, reg_democracy>30), rdrobust(y=deficit_perc_avg , x=businesswinmargin, bwselect="mserd",
                                                     covs=cbind(factor(unit_type),
                                                                deficit_perc_year0,
                                                                totalexpend_log_year0
                                                     ), all=TRUE, kernel='uni',p=1,cluster=region))

est3<-with(subset(vrn_b, reg_democracy<=30), rdrobust(y=competitive_all , x=businesswinmargin, bwselect="mserd",
                                                      covs=cbind(factor(unit_type),
                                                                 totalexpend_log_year0
                                                      ), all=TRUE, kernel='uni',p=1,cluster=region))

est4<-with(subset(vrn_b, reg_democracy>30), rdrobust(y=competitive_all , x=businesswinmargin, bwselect="mserd",
                                                     covs=cbind(factor(unit_type),
                                                                totalexpend_log_year0
                                                     ), all=TRUE, kernel='uni',p=1,cluster=region))

est5<-with(subset(vrn_b, reg_democracy<=30), rdrobust(y=competitive_construction , x=businesswinmargin, bwselect="mserd",
                                                      covs=cbind(factor(unit_type),
                                                                 totalexpend_log_year0
                                                      ), all=TRUE, kernel='uni',p=1,cluster=region))

est6<-with(subset(vrn_b, reg_democracy>30), rdrobust(y=competitive_construction , x=businesswinmargin, bwselect="mserd",
                                                     covs=cbind(factor(unit_type),
                                                                totalexpend_log_year0
                                                     ), all=TRUE, kernel='uni',p=1,cluster=region))



###########################################
########     TABLE   E5      ################
###########################################

#### Panel A

est1<-with(subset(vrn_b, party_alignment==0), rdrobust(y=natecon_avg , x=businesswinmargin, bwselect="mserd",
                                                       covs=cbind(factor(unit_type),
                                                                  natecon_year0,
                                                                  totalexpend_log_year0
                                                       ), all=TRUE, kernel='uni',p=1,cluster=region))

est2<-with(subset(vrn_b, party_alignment==1), rdrobust(y=natecon_avg , x=businesswinmargin, bwselect="mserd",
                                                       covs=cbind(factor(unit_type),
                                                                  natecon_year0,
                                                                  totalexpend_log_year0
                                                       ), all=TRUE, kernel='uni',p=1,cluster=region))

est3<-with(subset(vrn_b, party_alignment==0), rdrobust(y=education_avg , x=businesswinmargin, bwselect="mserd",
                                                       covs=cbind(factor(unit_type),
                                                                  education_year0,
                                                                  totalexpend_log_year0
                                                       ), all=TRUE, kernel='uni',p=1,cluster=region))

est4<-with(subset(vrn_b, party_alignment==1), rdrobust(y=education_avg , x=businesswinmargin, bwselect="mserd",
                                                       covs=cbind(factor(unit_type),
                                                                  education_year0,
                                                                  totalexpend_log_year0
                                                       ), all=TRUE, kernel='uni',p=1,cluster=region))

est5<-with(subset(vrn_b, party_alignment==0), rdrobust(y=health_avg , x=businesswinmargin, bwselect="mserd",
                                                       covs=cbind(factor(unit_type),
                                                                  health_year0,
                                                                  totalexpend_log_year0
                                                       ), all=TRUE, kernel='uni',p=1,cluster=region))

est6<-with(subset(vrn_b, party_alignment==1), rdrobust(y=health_avg , x=businesswinmargin, bwselect="mserd",
                                                       covs=cbind(factor(unit_type),
                                                                  health_year0,
                                                                  totalexpend_log_year0
                                                       ), all=TRUE, kernel='uni',p=1,cluster=region))

est7<-with(subset(vrn_b, party_alignment==0), rdrobust(y=totalexpend_avg_log , x=businesswinmargin, bwselect="mserd",
                                                       covs=cbind(factor(unit_type),
                                                                  totalexpend_log_year0
                                                       ), all=TRUE, kernel='uni',p=1,cluster=region))

est8<-with(subset(vrn_b, party_alignment==1), rdrobust(y=totalexpend_avg_log , x=businesswinmargin, bwselect="mserd",
                                                       covs=cbind(factor(unit_type),
                                                                  totalexpend_log_year0
                                                       ), all=TRUE, kernel='uni',p=1,cluster=region))


#### Panel B

est1<-with(subset(vrn_b, party_alignment==0), rdrobust(y=deficit_perc_avg , x=businesswinmargin, bwselect="mserd",
                                                       covs=cbind(factor(unit_type),
                                                                  deficit_perc_year0,
                                                                  totalexpend_log_year0
                                                       ), all=TRUE, kernel='uni',p=1,cluster=region))

est2<-with(subset(vrn_b, party_alignment==1), rdrobust(y=deficit_perc_avg , x=businesswinmargin, bwselect="mserd",
                                                       covs=cbind(factor(unit_type),
                                                                  deficit_perc_year0,
                                                                  totalexpend_log_year0
                                                       ), all=TRUE, kernel='uni',p=1,cluster=region))

est3<-with(subset(vrn_b, party_alignment==0), rdrobust(y=competitive_all , x=businesswinmargin, bwselect="mserd",
                                                       covs=cbind(factor(unit_type),
                                                                  totalexpend_log_year0
                                                       ), all=TRUE, kernel='uni',p=1,cluster=region))

est4<-with(subset(vrn_b, party_alignment==1), rdrobust(y=competitive_all , x=businesswinmargin, bwselect="mserd",
                                                       covs=cbind(factor(unit_type),
                                                                  totalexpend_log_year0
                                                       ), all=TRUE, kernel='uni',p=1,cluster=region))

est5<-with(subset(vrn_b, party_alignment==0), rdrobust(y=competitive_construction , x=businesswinmargin, bwselect="mserd",
                                                       covs=cbind(factor(unit_type),
                                                                  totalexpend_log_year0
                                                       ), all=TRUE, kernel='uni',p=1,cluster=region))

est6<-with(subset(vrn_b, party_alignment==1), rdrobust(y=competitive_construction , x=businesswinmargin, bwselect="mserd",
                                                       covs=cbind(factor(unit_type),
                                                                  totalexpend_log_year0
                                                       ), all=TRUE, kernel='uni',p=1,cluster=region))


###########################################
########     TABLE   E6      ################
###########################################

#### Panel A

est1<-with(subset(vrn_b, ur==0), rdrobust(y=natecon_avg , x=businesswinmargin, bwselect="mserd",
                                          covs=cbind(factor(unit_type),
                                                     natecon_year0,
                                                     totalexpend_log_year0
                                          ), all=TRUE, kernel='uni',p=1,cluster=region))

est2<-with(subset(vrn_b, ur==1), rdrobust(y=natecon_avg , x=businesswinmargin, bwselect="mserd",
                                          covs=cbind(factor(unit_type),
                                                     natecon_year0,
                                                     totalexpend_log_year0
                                          ), all=TRUE, kernel='uni',p=1,cluster=region))

est3<-with(subset(vrn_b, ur==0), rdrobust(y=education_avg , x=businesswinmargin, bwselect="mserd",
                                          covs=cbind(factor(unit_type),
                                                     education_year0,
                                                     totalexpend_log_year0
                                          ), all=TRUE, kernel='uni',p=1,cluster=region))

est4<-with(subset(vrn_b, ur==1), rdrobust(y=education_avg , x=businesswinmargin, bwselect="mserd",
                                          covs=cbind(factor(unit_type),
                                                     education_year0,
                                                     totalexpend_log_year0
                                          ), all=TRUE, kernel='uni',p=1,cluster=region))

est5<-with(subset(vrn_b, ur==0), rdrobust(y=health_avg , x=businesswinmargin, bwselect="mserd",
                                          covs=cbind(factor(unit_type),
                                                     health_year0,
                                                     totalexpend_log_year0
                                          ), all=TRUE, kernel='uni',p=1,cluster=region))

est6<-with(subset(vrn_b, ur==1), rdrobust(y=health_avg , x=businesswinmargin, bwselect="mserd",
                                          covs=cbind(factor(unit_type),
                                                     health_year0,
                                                     totalexpend_log_year0
                                          ), all=TRUE, kernel='uni',p=1,cluster=region))

est7<-with(subset(vrn_b, ur==0), rdrobust(y=totalexpend_avg_log , x=businesswinmargin, bwselect="mserd",
                                          covs=cbind(factor(unit_type),
                                                     totalexpend_log_year0
                                          ), all=TRUE, kernel='uni',p=1,cluster=region))

est8<-with(subset(vrn_b, ur==1), rdrobust(y=totalexpend_avg_log , x=businesswinmargin, bwselect="mserd",
                                          covs=cbind(factor(unit_type),
                                                     totalexpend_log_year0
                                          ), all=TRUE, kernel='uni',p=1,cluster=region))

### Panel B

est1<-with(subset(vrn_b, ur==0), rdrobust(y=deficit_perc_avg , x=businesswinmargin, bwselect="mserd",
                                          covs=cbind(factor(unit_type),
                                                     deficit_perc_year0,
                                                     totalexpend_log_year0
                                          ), all=TRUE, kernel='uni',p=1,cluster=region))

est2<-with(subset(vrn_b, ur==1), rdrobust(y=deficit_perc_avg , x=businesswinmargin, bwselect="mserd",
                                          covs=cbind(factor(unit_type),
                                                     deficit_perc_year0,
                                                     totalexpend_log_year0
                                          ), all=TRUE, kernel='uni',p=1,cluster=region))

est3<-with(subset(vrn_b, ur==0), rdrobust(y=competitive_all , x=businesswinmargin, bwselect="mserd",
                                          covs=cbind(factor(unit_type),
                                                     totalexpend_log_year0
                                          ), all=TRUE, kernel='uni',p=1,cluster=region))

est4<-with(subset(vrn_b, ur==1), rdrobust(y=competitive_all , x=businesswinmargin, bwselect="mserd",
                                          covs=cbind(factor(unit_type),
                                                     totalexpend_log_year0
                                          ), all=TRUE, kernel='uni',p=1,cluster=region))

est5<-with(subset(vrn_b, ur==0), rdrobust(y=competitive_construction , x=businesswinmargin, bwselect="mserd",
                                          covs=cbind(factor(unit_type),
                                                     totalexpend_log_year0
                                          ), all=TRUE, kernel='uni',p=1,cluster=region))

est6<-with(subset(vrn_b, ur==1), rdrobust(y=competitive_construction , x=businesswinmargin, bwselect="mserd",
                                          covs=cbind(factor(unit_type),
                                                     totalexpend_log_year0
                                          ), all=TRUE, kernel='uni',p=1,cluster=region))

###########################################
########     FIGURE E1      ################
###########################################

### Each set of regressions corresponds to the five annual regressions

natecon1<-with(vrn_b, rdrobust(y=natecon_year1 , x=businesswinmargin, bwselect="mserd",
                               covs=cbind(factor(unit_type), natecon_year0,totalexpend_log_year0), all=TRUE, kernel='uni',p=1,cluster=year))

natecon2<-with(vrn_b, rdrobust(y=natecon_year2 , x=businesswinmargin, bwselect="mserd",
                               covs=cbind(factor(unit_type), natecon_year0,totalexpend_log_year0), all=TRUE, kernel='uni',p=1,cluster=year))

natecon3<-with(vrn_b, rdrobust(y=natecon_year3 , x=businesswinmargin, bwselect="mserd",
                               covs=cbind(factor(unit_type), natecon_year0,totalexpend_log_year0), all=TRUE, kernel='uni',p=1,cluster=year))

natecon4<-with(vrn_b, rdrobust(y=natecon_year4 , x=businesswinmargin, bwselect="mserd",
                               covs=cbind(factor(unit_type), natecon_year0,totalexpend_log_year0), all=TRUE, kernel='uni',p=1,cluster=year))

natecon5<-with(vrn_b, rdrobust(y=natecon_year5 , x=businesswinmargin, bwselect="mserd",
                               covs=cbind(factor(unit_type), natecon_year0,totalexpend_log_year0), all=TRUE, kernel='uni',p=1,cluster=year))

### Health

health1<-with(vrn_b, rdrobust(y=health_year1 , x=businesswinmargin, bwselect="mserd",
                              covs=cbind(factor(unit_type), health_year0,totalexpend_log_year0), all=TRUE, kernel='uni',p=1,cluster=year))

health2<-with(vrn_b, rdrobust(y=health_year2 , x=businesswinmargin, bwselect="mserd",
                              covs=cbind(factor(unit_type), health_year0,totalexpend_log_year0), all=TRUE, kernel='uni',p=1,cluster=year))

health3<-with(vrn_b, rdrobust(y=health_year3 , x=businesswinmargin, bwselect="mserd",
                              covs=cbind(factor(unit_type), health_year0,totalexpend_log_year0), all=TRUE, kernel='uni',p=1,cluster=year))

health4<-with(vrn_b, rdrobust(y=health_year4 , x=businesswinmargin, bwselect="mserd",
                              covs=cbind(factor(unit_type), health_year0,totalexpend_log_year0), all=TRUE, kernel='uni',p=1,cluster=year))

health5<-with(vrn_b, rdrobust(y=health_year5 , x=businesswinmargin, bwselect="mserd",
                              covs=cbind(factor(unit_type), health_year0,totalexpend_log_year0), all=TRUE, kernel='uni',p=1,cluster=year))

#### Education

education1<-with(vrn_b, rdrobust(y=education_year1 , x=businesswinmargin, bwselect="mserd",
                                 covs=cbind(factor(unit_type), education_year0,totalexpend_log_year0), all=TRUE, kernel='uni',p=1,cluster=year))

education2<-with(vrn_b, rdrobust(y=education_year2 , x=businesswinmargin, bwselect="mserd",
                                 covs=cbind(factor(unit_type), education_year0,totalexpend_log_year0), all=TRUE, kernel='uni',p=1,cluster=year))

education3<-with(vrn_b, rdrobust(y=education_year3 , x=businesswinmargin, bwselect="mserd",
                                 covs=cbind(factor(unit_type), education_year0,totalexpend_log_year0), all=TRUE, kernel='uni',p=1,cluster=year))

education4<-with(vrn_b, rdrobust(y=education_year4 , x=businesswinmargin, bwselect="mserd",
                                 covs=cbind(factor(unit_type), education_year0,totalexpend_log_year0), all=TRUE, kernel='uni',p=1,cluster=year))

education5<-with(vrn_b, rdrobust(y=education_year5 , x=businesswinmargin, bwselect="mserd",
                                 covs=cbind(factor(unit_type), education_year0,totalexpend_log_year0), all=TRUE, kernel='uni',p=1,cluster=year))

#### Deficit 

deficit_perc1<-with(vrn_b, rdrobust(y=deficit_perc_year1 , x=businesswinmargin, bwselect="mserd",
                                    covs=cbind(factor(unit_type), deficit_perc_year0,totalexpend_log_year0), all=TRUE, kernel='uni',p=1,cluster=year))

deficit_perc2<-with(vrn_b, rdrobust(y=deficit_perc_year2 , x=businesswinmargin, bwselect="mserd",
                                    covs=cbind(factor(unit_type), deficit_perc_year0,totalexpend_log_year0), all=TRUE, kernel='uni',p=1,cluster=year))

deficit_perc3<-with(vrn_b, rdrobust(y=deficit_perc_year3 , x=businesswinmargin, bwselect="mserd",
                                    covs=cbind(factor(unit_type), deficit_perc_year0,totalexpend_log_year0), all=TRUE, kernel='uni',p=1,cluster=year))

deficit_perc4<-with(vrn_b, rdrobust(y=deficit_perc_year4 , x=businesswinmargin, bwselect="mserd",
                                    covs=cbind(factor(unit_type), deficit_perc_year0,totalexpend_log_year0), all=TRUE, kernel='uni',p=1,cluster=year))

deficit_perc5<-with(vrn_b, rdrobust(y=deficit_perc_year5 , x=businesswinmargin, bwselect="mserd",
                                    covs=cbind(factor(unit_type), deficit_perc_year0,totalexpend_log_year0), all=TRUE, kernel='uni',p=1,cluster=year))

totalexpend_log1<-with(vrn_b, rdrobust(y=totalexpend_log_year1 , x=businesswinmargin, bwselect="mserd",
                                       covs=cbind(factor(unit_type), totalexpend_log_year0), all=TRUE, kernel='uni',p=1,cluster=year))

totalexpend_log2<-with(vrn_b, rdrobust(y=totalexpend_log_year2 , x=businesswinmargin, bwselect="mserd",
                                       covs=cbind(factor(unit_type), totalexpend_log_year0), all=TRUE, kernel='uni',p=1,cluster=year))

totalexpend_log3<-with(vrn_b, rdrobust(y=totalexpend_log_year3 , x=businesswinmargin, bwselect="mserd",
                                       covs=cbind(factor(unit_type), totalexpend_log_year0), all=TRUE, kernel='uni',p=1,cluster=year))

totalexpend_log4<-with(vrn_b, rdrobust(y=totalexpend_log_year4 , x=businesswinmargin, bwselect="mserd",
                                       covs=cbind(factor(unit_type), totalexpend_log_year0), all=TRUE, kernel='uni',p=1,cluster=year))

totalexpend_log5<-with(vrn_b, rdrobust(y=totalexpend_log_year5 , x=businesswinmargin, bwselect="mserd",
                                       covs=cbind(factor(unit_type), totalexpend_log_year0), all=TRUE, kernel='uni',p=1,cluster=year))


